Soluciones Practico 1

Ejercicio 1

Seteo el ambiente

import numpy as np
import matplotlib.pyplot as plt
plt.style.use('ggplot')

Parte a

Escribir un programa mediante el cual simular y visualizar la evolución del capital del jugador A

def apuesta(p:float):
  if np.random.random() < p:
    return 1
  else:
    return -1

def juego(p:float, S:int, X0:int):
  trayectoria = [X0]
  X = X0
  while (X > 0) and (X < S):
    X += apuesta(p)
    trayectoria.append(X)
  return trayectoria

p = 0.5
S = 10
X0 = 5
trayectoria = juego(p, S, X0)

print('Trayectoria: ', trayectoria)

Trayectoria:  [5, 4, 3, 4, 3, 2, 1, 2, 1, 2, 3, 4, 3, 2, 3, 2, 1, 0]

plt.figure(figsize=(8,5))
plt.plot(trayectoria, marker='o')
plt.title('Juego de la ruina - Una trayectoria')
plt.ylim(-0.5, S+0.5)
plt.axhline(y=S, linestyle='--')
plt.axhline(y=0, linestyle='--')
plt.xlabel(r'Turno $n$')
plt.ylabel(r'Capital $X_n$')
plt.show()

todas_trayectorias = []
for _ in range(5):
  todas_trayectorias.append(juego(p, S, X0))

plt.figure(figsize=(8,5))
for trayectoria in todas_trayectorias:
  plt.plot(trayectoria, marker='o')
plt.title('Juego de la ruina - Varias trayectorias') # Added title for clarity
plt.ylim(-0.5, S+0.5)
plt.axhline(y=S, linestyle='--')
plt.axhline(y=0, linestyle='--')
plt.xlabel(r'Turno $n$')
plt.ylabel(r'Capital $X_n$')
plt.show()

Parte b

Desarrollar el programa anterior de manera de estimar, mediante la simulación de un número adecuado de trayectorias del capital del jugador A, la probabilidad de ruina.

def probabilidad_ruina(p:float, S:int, X0:int, N:int):
  ruinas = 0
  for _ in range(N):
    trayectoria = juego(p, S, X0)
    if trayectoria[-1] == 0:
      ruinas += 1
  return ruinas/N

p = 0.7
S = 10
X0 = 5
N = 10000
print(f'Probabilidad de ruina: {probabilidad_ruina(p, S, X0, N)}')

Probabilidad de ruina: 0.0141

fS = []
for X0 in range(S+1):
  fS.append(probabilidad_ruina(p, S, X0, N))

plt.figure(figsize=(8,5))
plt.plot(fS, marker='o')
plt.title('Probabilidad de ruina')
plt.xlabel(r'Capital inicial $X_0$')
plt.show()

(Extra) Propiedad de Markov

def probabilidad_ruina_Markov(p:float, S:int, X0:int, X1:int, N:int):
  ruinas = 0
  condicion = 0
  for _ in range(N):
    trayectoria = juego(p, S, X0)
    if trayectoria[1] == X1:
      condicion += 1
      if trayectoria[-1] == 0:
        ruinas += 1
  return ruinas/condicion

p = 0.5
S = 10
X0 = 5
N = 10000
p0 = probabilidad_ruina(p, S, X0+1, N)
p1 = probabilidad_ruina_Markov(p, S, X0, X0+1, N)
print(f'p0: {p0}\np1: {p1}')

p0: 0.3962
p1: 0.39119222376512597

Ejercicio 2

Parte a

Estimación de la duración esperada del juego

def simular_duracion_juego(p: float, S: int, X0: int) -> int:
    trayectoria = [X0]
    X = X0
    turnos = 0
    while (X > 0) and (X < S):
        X += apuesta(p) ## funcion del ejercicio 1
        turnos += 1
    return turnos

def estimar_duracion_esperada(p: float, S: int, X0: int, N: int) -> float:
    duraciones = []
    for _ in range(N):
        duraciones.append(simular_duracion_juego(p, S, X0))
    return np.mean(duraciones)

(i) Duración esperada vs. capital inicial \(X_0 = k\)

p_ej2_i = 0.5
S_ej2_i = 10
N_ej2_i = 5000 

duracion_vs_X0 = []
for k_val in range(S_ej2_i + 1):
    duracion_vs_X0.append(estimar_duracion_esperada(p_ej2_i, S_ej2_i, k_val, N_ej2_i))

print("Duración esperada vs. Capital inicial:", [float(d) for d in duracion_vs_X0])

Duración esperada vs. Capital inicial: [0.0, 8.5572, 16.3292, 20.8324, 23.7448, 25.1652, 24.4272, 20.8436, 15.5248, 9.3192, 0.0]

plt.figure(figsize=(8,5))
plt.plot(range(S_ej2_i + 1), duracion_vs_X0, marker='o')
plt.title(f'Duración Esperada vs. Capital Inicial (S={S_ej2_i}, p={p_ej2_i})')
plt.xlabel(r'Capital inicial $X_0$')
plt.ylabel('Duración Esperada del Juego')
plt.xticks(range(S_ej2_i + 1))
plt.grid(True)
plt.show()

(ii) Duración esperada vs. probabilidad \(p\)

X0_ej2_ii = 5
S_ej2_ii = 10
N_ej2_ii = 5000 

probabilidades = np.linspace(0.1, 0.9, 9) 
duracion_vs_p = []
for p_val in probabilidades:
    duracion_vs_p.append(estimar_duracion_esperada(p_val, S_ej2_ii, X0_ej2_ii, N_ej2_ii))

print("Duración esperada vs. Probabilidad p:", [float(d) for d in duracion_vs_p])

Duración esperada vs. Probabilidad p: [6.2168, 8.254, 12.3132, 19.2708, 25.3448, 19.6704, 12.3836, 8.302, 6.2856]

plt.figure(figsize=(8,5))
plt.plot(probabilidades, duracion_vs_p, marker='o')
plt.title(f'Duración Esperada vs. Probabilidad p (S={S_ej2_ii}, X0={X0_ej2_ii})')
plt.xlabel(r'Probabilidad $p$ de ganar')
plt.ylabel('Duración Esperada del Juego')
plt.grid(True)
plt.show()

Parte b

Visualizar las estimaciones de la duración esperada del juego y compararlas con las correspondientes expresiones analíticas.

Funciones analíticas

def prob_ruina_analitica(p: float, S: int, k: int) -> float:
    q = 1 - p
    if p == 0.5:
        return (S - k) / S
    else:
        ratio = q / p
        return (ratio**k - ratio**S) / (1 - ratio**S)

def duracion_esperada_analitica(p: float, S: int, k: int) -> float:
    q = 1 - p
    if p == 0.5:
        return k * (S - k)
    else:
        prob_ganar_analitica = 1 - prob_ruina_analitica(p, S, k)
        if p != q:
            return (k / (q - p)) - (S / (q - p)) * ((1 - (q / p)**k) / (1 - (q / p)**S))
        else:
            return k * (S - k)

(i) Duración esperada vs. capital inicial \(X_0 = k\) (Comparación)

duracion_analitica_X0 = []
for k_val in range(S_ej2_i + 1):
    duracion_analitica_X0.append(duracion_esperada_analitica(p_ej2_i, S_ej2_i, k_val))

print("Duración analítica vs. Capital inicial:", [float(d) for d in duracion_analitica_X0])

Duración analítica vs. Capital inicial: [0.0, 9.0, 16.0, 21.0, 24.0, 25.0, 24.0, 21.0, 16.0, 9.0, 0.0]

plt.figure(figsize=(10,6)) 
plt.plot(range(S_ej2_i + 1), duracion_vs_X0, marker='o', label='Simulado', linestyle='--')
plt.plot(range(S_ej2_i + 1), duracion_analitica_X0, marker='x', label='Analítico', linestyle='-')
plt.title(f'Comparación Duración Esperada vs. Capital Inicial (S={S_ej2_i}, p={p_ej2_i})')
plt.xlabel(r'Capital inicial $X_0$')
plt.ylabel('Duración Esperada del Juego')
plt.xticks(range(S_ej2_i + 1))
plt.grid(True)
plt.legend()
plt.show()

(ii) Duración esperada vs. probabilidad \(p\) (Comparación)

duracion_analitica_p = []
for p_val in probabilidades:
    duracion_analitica_p.append(duracion_esperada_analitica(p_val, S_ej2_ii, X0_ej2_ii))

print("Duración analítica vs. Probabilidad p:", [float(d) for d in duracion_analitica_p])

Duración analítica vs. Probabilidad p: [6.249788314987299, 8.317073170731707, 12.14369501466276, 19.181818181818183, 25.0, 19.18181818181818, 12.143695014662754, 8.317073170731707, 6.249788314987299]

plt.figure(figsize=(10,6))
plt.plot(probabilidades, duracion_vs_p, marker='o', label='Simulado', linestyle='--')
plt.plot(probabilidades, duracion_analitica_p, marker='x', label='Analítico', linestyle='-')
plt.title(f'Comparación Duración Esperada vs. Probabilidad p (S={S_ej2_ii}, X0={X0_ej2_ii})')
plt.xlabel(r'Probabilidad $p$ de ganar')
plt.ylabel('Duración Esperada del Juego')
plt.grid(True)
plt.legend()
plt.show()

Ejercicio 3

Consideramos un juego de ruina del jugador con posibilidad de empate (también llamado lazy random walk): en el tiempo \(n\) la ganancia \(X_n\) del Jugador A puede aumentar en una unidad con probabilidad \(r \in (0; 1/2]\), disminuir en una unidad con probabilidad \(r\), o permanecer igual con probabilidad \(1 - 2r\). Definimos

\[ f (k) := P(RA | X_0 = k) \]

como la probabilidad de ruina del Jugador A, y

\[ h(k) := E[T_{0;S} | X_0 = k] \]

como la esperanza de la duración del juego \(T_{0;S}\) comenzando en \(X_0 = k\), para \(k = 0; 1; \dots ; S\).

a) Ecuación en Diferencias para la Probabilidad de Ruina \(u_k\)

Para encontrar la ecuación que satisface \(u_k\), la probabilidad de que el Jugador A se arruine comenzando con una fortuna k, utilizamos un análisis de primer paso. Consideramos el resultado del primer movimiento desde el estado k:

• Gana: Con probabilidad \(p\), la fortuna se convierte en \(k+1\). La probabilidad de ruina desde este nuevo estado es \(u_{k+1}\).
• Pierde: Con probabilidad \(q\), la fortuna se convierte en \(k-1\). La probabilidad de ruina desde este nuevo estado es \(u_{k-1}\).
• Empata: Con probabilidad \(r\), la fortuna permanece en \(k\). La probabilidad de ruina sigue siendo \(u_k\).

Por la ley de probabilidad total, podemos escribir: \(u_k = p \cdot u_{k+1} + q \cdot u_{k-1} + r \cdot u_k\)

Para simplificar, agrupamos los términos de \(u_k\): \(u_k (1-r) = p \cdot u_{k+1} + q \cdot u_{k-1}\)

Dado que \(p+q+r=1\), tenemos que \(1-r=p+q\). Sustituyendo: \((p+q)u_k = p \cdot u_{k+1} + q \cdot u_{k-1}\)

Reordenando obtenemos la ecuación en diferencias homogénea:

\(pu_{k+1} - (p+q)u_k + qu_{k-1} = 0, \text{para } 1 \leq k \leq N-1\) Las condiciones de contorno son los estados absorbentes:

• Si el jugador comienza en \(k=0\), ya está en la ruina: \(u_0=1\).
• Si el jugador comienza en \(k=N\), ha ganado todo, por lo que no puede arruinarse: \(u_N=0\).

\(u_0=1, u_N=0\)

b) Solución de la Ecuación Homogénea

Resolvemos la ecuación mediante su polinomio característico: \(p\lambda^2 - (p+q)\lambda + q = 0\)

Este polinomio se puede factorizar como \((p\lambda-q)(\lambda-1)=0\). Las raíces son: \(\lambda_1=1 \text{ y } \lambda_2 = \frac{q}{p}\).

Tenemos dos casos:

Caso 1: \(p \neq q\) Las raíces son distintas. La solución general es de la forma: \(u_k = A(\lambda_1)^k + B(\lambda_2)^k = A+B\left(\frac{q}{p}\right)^k\)

Aplicamos las condiciones de contorno:

\(u_0 = 1 \implies A+B=1\)

\(u_N = 0 \implies A+B\left(\frac{q}{p}\right)^N = 0\)

Resolviendo el sistema, obtenemos \(B = \frac{1}{1-(q/p)^N}\) y \(A = -B(q/p)^N\). La solución es:

\(u_k = \frac{(q/p)^k - (q/p)^N}{1-(q/p)^N}\)

Caso 2: \(p=q\) Las raíces son iguales (\(\lambda_1=\lambda_2=1\)). La solución general es: \(u_k = A(1)^k + Bk(1)^k = A+Bk\)

Aplicamos las condiciones de contorno:

\(u_0 = 1 \implies A=1\)

\(u_N = 0 \implies A+BN=0 \implies 1+BN=0 \implies B = -\frac{1}{N}\)

La solución es:

\(u_k = 1-\frac{k}{N} = \frac{N-k}{N}\)

Compatibilidad con la intuición: Sí, la solución es compatible. Notablemente, la probabilidad de ruina no depende de \(r\) (la probabilidad de empate). Un empate simplemente retrasa el juego, pero no altera el sesgo fundamental (\(p\) vs \(q\)) que determina el resultado final. La probabilidad de ruina es idéntica a la del juego de ruina del jugador estándar sin empates.

c) Ecuación en Diferencias para la Duración Media \(d_k\)

Nuevamente, usamos un análisis de primer paso para \(d_k\), la duración media del juego comenzando en \(k\). Cada paso toma 1 unidad de tiempo. \(d_k = 1 + p \cdot d_{k+1} + q \cdot d_{k-1} + r \cdot d_k\)

Agrupando los términos de \(d_k\) y usando \(1-r=p+q\): \((1-r)d_k = 1 + p \cdot d_{k+1} + q \cdot d_{k-1}\) \((p+q)d_k = 1 + p \cdot d_{k+1} + q \cdot d_{k-1}\)

Reordenando, obtenemos la ecuación en diferencias no homogénea:

\(pd_{k+1} - (p+q)d_k + qd_{k-1} = -1, \text{para } 1 \leq k \leq N-1\) Las condiciones de contorno son: Si el juego comienza en un estado absorbente, la duración es 0.

\(d_0=0, d_N=0\)

d) Solución Particular de la Ecuación para \(d_k\)

Buscamos una solución particular \(d_k^{(p)}\) para la ecuación no homogénea.

Caso 1: \(p \neq q\) Intentamos con una solución de la forma \(d_k^{(p)}=Ck\). \(pC(k+1)-(p+q)Ck+qC(k-1)=-1\) \(pCk+pC-pCk-qCk+qCk-qC=-1\) \(C(p-q)=-1 \implies C = \frac{1}{q-p}\)

Una solución particular es:

\(d_k^{(p)} = \frac{k}{q-p}\)

Caso 2: \(p=q\) El método anterior falla porque \(Ck\) es parte de la solución homogénea. Intentamos con \(d_k^{(p)} = Ck^2\). La ecuación es \(pd_{k+1}-2pd_k+pd_{k-1}=-1\). \(pC(k+1)^2 - 2pCk^2 + pC(k-1)^2 = -1\) \(pC(k^2+2k+1 - 2k^2 + k^2-2k+1) = -1\) \(pC(2)=-1 \implies C = -\frac{1}{2p}\)

Una solución particular es:

\(d_k^{(p)} = -\frac{k^2}{2p}\)

e) Solución Completa para \(d_k\)

La solución general es la suma de la solución de la homogénea (\(d_k^{(h)}\)) y una particular (\(d_k^{(p)}\)).

Caso 1: \(p \neq q\) \(d_k = d_k^{(h)} + d_k^{(p)} = A+B\left(\frac{q}{p}\right)^k + \frac{k}{q-p}\)

Aplicamos las condiciones de contorno (\(d_0=0, d_N=0\)):

\(d_0=0 \implies A+B=0 \implies A=-B\)

\(d_N=0 \implies A+B\left(\frac{q}{p}\right)^N + \frac{N}{q-p} = 0\)

Sustituyendo \(A=-B\): \(B\left(\left(\frac{q}{p}\right)^N - 1\right) = -\frac{N}{q-p} \implies B = \frac{N}{(p-q)((q/p)^N-1)}\)

La solución completa es:

\(d_k = \frac{k}{q-p} + \frac{N}{(p-q)(1-(q/p)^N)}\left(\left(\frac{q}{p}\right)^k - 1\right)\)

Caso 2: \(p=q\) \(d_k = d_k^{(h)} + d_k^{(p)} = A+Bk-\frac{k^2}{2p}\)

Aplicamos las condiciones de contorno (\(d_0=0, d_N=0\)):

\(d_0=0 \implies A=0\)

\(d_N=0 \implies BN - \frac{N^2}{2p} = 0 \implies B = \frac{N}{2p}\)

La solución completa es:

\(d_k = \frac{Nk}{2p} - \frac{k^2}{2p} = \frac{k(N-k)}{2p}\)

f) Comportamiento de la Duración Media \(d_k\) cuando \(p=q\)

Para el caso simétrico (\(p=q\)), la duración media es:

\(d_k = \frac{k(N-k)}{2p}\)

Recordando que \(p+q+r=1\) y \(p=q\), tenemos \(2p+r=1 \implies 2p=1-r\). La fórmula se puede reescribir como:

\(d_k = \frac{k(N-k)}{1-r}\)

Compatibilidad con la intuición: Sí, esta solución es muy intuitiva por dos razones principales:

• Forma de la función: La función \(k(N-k)\) es una parábola invertida que es 0 en \(k=0\) y \(k=N\). Alcanza su máximo en \(k=N/2\). Esto significa que la duración esperada del juego es mayor cuando el jugador comienza a mitad de camino, lo más lejos posible de ambos estados absorbentes (ruina o victoria total).

• Dependencia de \(r\): La duración media \(d_k\) es inversamente proporcional a \(1-r\). A medida que la probabilidad de empate \(r\) se acerca a 1, el denominador \(1-r\) se acerca a 0, y la duración media \(d_k \to \infty\). Esto tiene perfecto sentido: si los jugadores casi siempre empatan, el juego se “congela” en su estado actual y tardará una cantidad de tiempo extremadamente larga en llegar a su fin.